a = [1.25 8.75 0.5 5.75 3 7.25];
b = [1.25 0.75 4.75 5 6.5 7.25];
d = [3 5 4 7 6 11];
x = [5 2];
y = [1 7];
for i = 1:6 
    for j = 1:2
        dij(i, j) = sqrt((x(j)-a(i)).^2 + (y(j)-b(i)).^2);
    end
end

% 技巧: 先确定乘积列向量的行数
f = [dij(:, 1); dij(:, 2)]; % 两列aa合并成一列, 作为目标函数的系数
Aeq =  [1 0 0 0 0 0  1 0 0 0 0 0;
        0 1 0 0 0 0  0 1 0 0 0 0;
        0 0 1 0 0 0  0 0 1 0 0 0;
        0 0 0 1 0 0  0 0 0 1 0 0;
        0 0 0 0 1 0  0 0 0 0 1 0;
        0 0 0 0 0 1  0 0 0 0 0 1];
beq = d.';
A = [1 1 1 1 1 1 0 0 0 0 0 0;
     0 0 0 0 0 0 1 1 1 1 1 1];
B = [20; 20]; % 列
LB = zeros(12, 1); % 列
UB = [];

[optx, fval] = linprog(f, A, B, Aeq, beq, LB, UB);

